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Abstract. Analytical representations are derived for two equations of state (EOSs) of neutron-star matter: FPS 
and SLy. Each of these EOSs is unified, that is, it describes the crust and the core of a neutron star using the 
same physical model. Two versions of the EOS parametrization are considered. In the first one, pressure and mass 
density are given as functions of the baryon density. In the second version, pressure, mass density, and baryon 
density are given as functions of the pseudo-enthalpy, which makes this representation particularly useful for 2-D 
calculations of stationary rotating configurations of neutron stars. 
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1. Introduction 

The equation of state (EOS) of dense matter is a crucial 
input for the neutron-star structure calculations. Under 
standard conditions neutron-star matter is strongly de- 
generate, and therefore the matter pressure is temperature 
independent; exceptions are the outermost (a few meters 
thick) envelopes, newly-born neutron stars, and the en- 
velopes of exploding X-ray bursters. At p <; 10 8 g cm -3 
the EOS is not affected by the magnetic field even as 
strong as 10 14 G and by the temperature T ^ 10 9 K. 
Therefore, except for a thin outer envelope (which for the 
bottom density of 10 s g cm~ 3 is ~ 100 m thick and con- 
tains only ~ 10 -9 M© of matter), the EOS of neutron-star 
matter has a one-parameter character. In order to deter- 
mine the neutron-star structure up to the maximum allow- 
able mass, M max , one has to know the EOS up to a few 
times 10 15 g cm" 3 . While the EOS of the neutron-star 
crust (p ^ (1 — 2) x 10 14 g cm -3 ) is rather well known 
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Ce.g.. lHaenselll200ll) . the EOS of the liquid core at 
5 x 10 14 g cm -3 , which is crucial for determining Mm,,, 
mains uncerta in (e.g., iHeiselberg fc Pandharipandel l2000; 
lHaensellEool . 

The EOS of the crust depends on the crust formation 
scenario. Two limiting cases are: cold catalyzed matter be- 
ing at the ground state at a fixed baryon density n, and the 
accreted crust formed via compression at T < 10 9 K from 
the thermonuclear ashes of the X-ray bursts in the outer 
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envelopes o f accreting neutron stars l|Haensel fc ZdunikL 
Il99(j l2003^> . The EOS of the liquid core does not depend 
on the formation scenario, and can slightly deviate from 
the ground state one only when there are deviations from 
the weak interaction equilibrium - e.g., when rapid matter 
flow (like in stellar pulsations) is involved. In the present 
paper we consider the standard case of the ground-state 
matter. 

A "unified EOS" is obtained in the many-body cal- 
culations based on a single effective nuclear Hamiltonian, 
and is valid in all regions of the neutron star interior. For 
unified EOSs the transitions between the outer crust and 
the inner crust, and between the inner crust and the core 
are obtained as a result of many-body calculation. Alas, 
up to now only a few models of unified EOSs have been 
constructed. All other EOSs consist of crust and core seg- 
ments obtained using different physical models. The crust- 
core interface has there no physical meaning and both 
segments are joined using an ad hoc matching procedure. 
Therefore, neutron-star models based on these EOSs con- 
tain a shell with an unphysical EOS. For such "matched 
EOSs" it is not possible to study phenomena which are 
sensitive to the position of the crust-core interface. 

In the pr esent paper we consider two unified E OSs, the 
FPS EOS o f iPandharipande fc Ravenhalll ^1989j) and the 
SLy EOS of lDouchin fc Haensell fcOOll) . 

EOSs are usually given in the form of tables. Therefore, 
in order to use them, one has to employ interpolation 
between the tabulated points. However, the interpola- 
tion procedure is not unique. This introduces ambiguities 
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in the calculated parameters of the neutron star mod- 
els. Moreover, interpolation should be done respecting 
exact thermodynamic relations. This turned out to be 
a particularly serious problem in the high-precision 2-D 
calculations of mode ls of rapidly rotating neutron stars 
ijNozawa et all Il998|) . In the 3-D calculations of the sta- 
tionary configurations of a close binary neutron-star sys- 
tem one needs derivatives of pressure with respect to the 
enthalpy, and tabulated non -polytropic EOSs are not use- 
ful in this respect (see, e.g.. lOourgoulhon et alllioflll) . 

In view of the deficiencies and ambiguities inherent 
in the tabulated EOSs, which are particularly serious in 
the case of matched EOSs, it is of great interest to derive 
analytical representations of the EOSs. They have two im- 
portant advantages over the tabulated ones. First, there 
is no ambiguity of interpolation and the derivatives can 
be precisely calculated. Second, these representations can 
be constructed fulfilling exactly the thermodynamic rela- 
tions. In this way, analytical EOSs can allow for a very 
high precision of neutron star structure calculation in the 
2-D and 3-D cases. 

In the present paper we derive analytical representa- 
tions of the unified FPS and SLy EOSs. In Sect. |2 we 
discuss general properties of the EOS used for the calcu- 
lations of non-rotating and rotating neutron-star models. 
Analytical representations of the EOSs for non-rotating 
stars are described in Sect. 03 while Sect. 0] is devoted to 
the case of rotating stars. Discussion and conclusion are 
presented in Sect. [51 

2. General properties of the one-parameter EOS 
of neutron-star matter 

A given EOS is usually presented in the form of a table 
containing a grid of calculated values of matter density 
(full energy density £ including the rest energies of the 
matter constituents divided by c 2 , i.e., p — £/c 2 ), baryon 
(number) density n, and pressure P. The EOS {pt, rii, Pi} 
(i = 1, . . . , N) is then interpolated between the tabulated 
points so as to get the one-parameter EOS in the form P = 
P(n), p = p(n). The interpolation involves some degree of 
indeterminacy (there are many ways of interpolating) and 
this itself implies some ambiguity as far as the calculated 
values of the neutron-star parameters (for example, the 
value of M max ) are concerned. 



2.1. Non-rotating configurations 

High precision determination of the baryon number A and 
gravitational mass M of an equilibrium configuration re- 
quires condition of thermodynamic consistency of func- 
tions p(n) and P(n) to be strictly fulfilled. The first law 
of thermodynamics in the T = limit implies relation 



can be also used in the integral forms, 



P(n) 



d 

dn 



(1) 



which p uts constraints on the inter polation procedure 
(see, e.g.. lrlaensel fc Pros zvnskilll982|) : the above relation 
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where p s and n s are the values of p and n at the neutron- 
star surface. In the present paper we put p s equal to 
density of 56 Fe at zero pressure and zero temperature, 
p s = 7.86 g cm~ 3 . In the outermost neutron-star layers, we 
fix the value of mass per nucleon as too = 1.66 x 10~ 24 g, 
so that n s = p s /m = 4.73494 x 10 24 cm" 3 . 

An example of an interpolation recipe which respects 
Eqs. CE])~© and at the same time yields highly smoot h 
functions P(n) and p(n) was presented bv lSwestvl ()l996|) . 
Only a thermodynamically consistent interpolation yields 
neutron-star models which strictly satisfy the relation con- 
necting the baryon chemical potential pb = c 2 dp/dnb = 
(pc 2 + P)/n and the metric function $(r) at a given cir- 
cumferential radius r. This strict relation stems from the 
equation for the metric function 



d$ 
dr 



pc* 



dP 



pc 2 + P dr 
and can be written as 

p b (r)e*M = Mb(i?)e* (fl) 



(4) 



(5) 



where $(i?) = y/l - 2GM/Rc 2 . Here, M and R are the 
total gravitational mass and circumferential stellar radius, 
respectively. Strictly speaking, if Eq. JSJ does not hold, the 
calculated configuration is not in hydrostatic equilibrium. 
Alas, this often happens for the neutron-star models calcu- 
lated using tabulated EOSs with logarithm ic interpolation 
between the tabulated points (see, e.g.. lHarrison et all 
ll96.4lBavm et al.ll97lUArnett & Bowersll977l) . Resulting 
inconsistencies may seem minor, but they may lead to se- 
rious problems if high precision of a simultaneous deter- 
mination of M and A is concerned; this is the case of the 
energy release due to a phase transition in the neutron- 
star core (see, e.g.. lHaensel fe ProszvnskiH l982). On the 
contrary, if Eq. JSJ is strictly fulfilled, then the constancy 
of pbe* and the accuracy of calculating stellar parame- 
ters are limited only by the numerical precision of the 
computer code. 

For a thermodynamically consistent EOS, Eq. JSJ im- 
plies the baryon density profile within a static neutron 
star, 



n(r) 



p(r) , P{r) 



Ps 



p s c z 



3 *(r)-*(fl) 



(6) 



where we have neglected P s = P(R, T) because the pres- 
sure within the atmosphere is small compared with p s c 2 . 
Therefore, what one needs to get n(r) in a static neutron 
star is just the surface baryon density, n s ; the {r^} column 
of the tabulated EOS turns out to be redundant. 
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2.2. Stationary rotating configurations 



Table 1. Parameters of the fit (TL4"|) 



Rotation breaks the spherical symmetry of the equilibrium 
configuration. Stationary configurations of a rigidly rotat- 
ing star of ideal liquid are solutions of 2-D axially sym- 
metric partial differential equations of hydrostatic equilib- 
rium in coordi nates r and 9 (rot ation of relativistic stars 
is reviewed by IStergioulasI 2003). As far as the EOS is 
concerned, it is suitable to parametrize it in terms of a 
dimcnsionless pseudo-enthalpy 



H(P) 



dP' 



p(P')c 2 + P' ' 
It can be rewritten in terms of the enthalpy per baryon 
pc 2 + P 



(7) 



h(P) = 



as 



H{P) = In 



h(P) 



(8) 



(9) 



where h s — p s c 2 /n s . The pseudo-enthalpy H vanishes at 
the stellar surface and increases monotonically towards 
the stellar center. It is a very useful variable in numerical 
calculations of the stationary configurations of rotating 
stars, because the definition of H allows one to easily write 
down the first integral of motion 

ff(r,0) + $(r,0) -ln 7 (r,0) = H(0, 0) + $(0, 0) , (10) 

where $ is the metric function and 7 is the Lorentz factor 
connecting the two measured values of a physical quantity: 
that measured in a local observer's reference frame at a 
point (r, 6, (f>) and that measured by the distant observer 
in an inertial reference frame at infinity, 



-1/2 



(11) 



where U is the fluid velocity in the azimuthal (<f>) direction, 
as measured by a local observer. Equation I|1U|) is the gen- 
eral relativistic version of the Bernoulli theorem for a sta- 
tionary ideal fluid flow in gravitational field. Consequently, 
the most useful parametrization of the EOS for rotating 
stars is p — p(H), P = P(H). Then, the baryon number 
density is given by an exact analytical formula 



P(H) P{H) 



-H 



(12) 



n{H) = n s 

p s p s C^ 

which is the rotating-star analogue of Eq. ©. 

3. Analytical representations of the EOS 

There are three qualitatively different domains of the in- 
terior of a neutron star, separated by phase transition 
points: the outer crust (consisting of the electrons and 
atomic nuclei), the inner crust (consisting of the electrons, 
nuclei, and dripped neutrons), and the core which contains 
the electrons, neutrons, protons, /i _ -mesons, and possi- 
bly 7T- and if-mesons, some hyperons, or quark matter. 
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The latter species are contained in the innermost stel- 
lar domain called the inner core. In addition, there can 
be density discontinuities at the interfaces between lay- 
ers containing different nuclei in the crust. In the fitting, 
we neglect these small discontinuities and approximate 
the EOS by fully analytical functions. However, the dif- 
ferent character of the EOS in the different domains is 
reflected by the complexity of the fit, which consists of 
several fractional-polynomial parts, matched together by 
virtue of the function 



fo(x) 



1 



e x + 1 



(13) 



We rely on a tabulated unified EOS (FPS 1 or SLy) at 
p > 5 x 10 10 g cm" 3 . At lower densities, 10 8 g cm" 3 < 
p < 5 x 10 1 g cm" 3 , the crus t al ma tter is described by 
the EOS of lHaensel fc Pichonl l|l994|) (HP94), based on 
experimental nuclear data, su pplemented by the E OS for 
cold catalyzed matter due to iBavm et al.l l)l97lh (BPS) 
at still lower density p Si 10 8 g cm~ 3 . The lowest-density 
parts of the tables at p < 10 5 g cm -3 have not been used 
in the fitting, because at such low density the EOS is no 
longer one-parametric, but depends also on temperature 
(see Fig.|TJ. 



3.1. Non-rotating stars 

For non-rotating configurations, we have parametrized 
the pressure as function of density. Let us denote 
£ = log(p/g cm" 3 ), C = log(P/dyn cm" 2 ). Then the 
parametrization reads 



ai +a 2 £ + a 3 £ 3 

C = r- ; /o(«5 £-a 6 

1 + 04 £ 

+ (a 7 + a 8 £) /o(ag(aio - 0) 

+ (an + ai 2 £) /o(ai3(ai4 - 0) 

+ (a i5 + ai 6 £) /o(ai7(ai8 - 0) ■ 



(14) 



The parameters for FPS and SLy EOSs are given in 
Table[T] The typical fit error of P is 1-2% (for £ ^ 5). The 
maximum error is determined by the jumps near the phase 
transitions in the tabulated EOSs, which are smoothed 



The FPS table has been kindly provided by N. Stergioulas. 
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Fig. 1. Neutron-star EOS for non-rotati o nal co nfigura- 
tions: BPS (triangles), lHaensel fc Pichonl <|l994 (HP94, 
stars), SLy (dots), OPAL at T = 10 6 , 10 7 , and 10 8 K 
(dashed lines), the fit IjMfl (solid line) and the fit modified 
at low p (dotted line). 



by the fit 114|) . For FPS, the maximum error is 3.6% at 
£ = 14.22 (crust-core interface). For SLy the maximum 
error is 2.9% at £ = 8.42 ( 62 Ni- 64 Ni interface in the HP94 
part of the table). 

Figure ^ shows logP against log/? for a tabulated 
EOS (symbols) and the corresponding fit (solid line). 
Triangles correspond to BPS, stars to HP94, and dots to 
SLy data. By construction, the fit is accurate at p 10 s 
g cm -3 . As stated above, at lower density the EOS be- 
comes temperature-dependent (for T values typical for 
neutron star envelopes). This is ill ustrated by the das hed 
lines, that show the OPAL EOS l|Rogers et all Il99f^ for 
T = 10 6 , 10 7 , and 10 s K. 2 However, a reasonable con- 
tinuation of the fit to lower densities can be constructed 
by a simple interpolation. For instance, the dotted line in 
Fig. ^ corresponds to P — Pa t + Po, where Pfit = 1CK 
is given by Eq. (|14|l (where we always assume £ > 0), 
and Po = 3.5 x 10 14 p approximates the OPAL EOS near 
p ~ p s at T = 10 7 K (here P is in dyn cm -2 and p in 
g cm -3 ). 

In Fig.|2we compare the FPS and SLy EOSs. Symbols 
on the upper panel show the data (triangles, stars, dots, 
and open circles for BPS, HP94, SLy, and FPS, respec- 
tively) and lines show the fit (solid for SLy and dot-dashed 
for FPS). In order to make the differences between the 
data and fits and between SLy and FPS EOSs visible, we 
plot the difference logP— 1.4 log p, where P is in dyn cm -2 
and p in g cm~ 3 . The lower panel shows the relative dif- 
ference between the tabulated and fitted EOSs (solid and 



2 The OPAL table for iron has been kindly provided by F. J. 
Rogers. 
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Fig. 2. Comparison of the data and fits for SLy and FPS 
EOSs for non-rotational configurations. Upper panel: rar- 
efied tabular data (symbols) and the fit l|14fl (lines); lower 
panel: relative difference between the data and fit. Filled 
dots and solid line: SLy; open circles and dot-dashed line: 
FPS (triangles and stars on the upper panel are BPS and 
HP94 data at p < 5 x 10 10 g cm -3 ). 

Table 2. Parameters of the fits (JTSJ and JTJJ 



i 


p s (FPS) 


Pi (SLy) 


gi(FPS) 


<?i(SLy) 


1 


0.320 


0.423 


0.608 


0.183 
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2.17 


2.42 


2.41 


1.26 
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0.173 


0.031 


2.39 


6.88 
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3.01 


0.78 


3.581 


3.612 


5 


0.540 


0.238 


1.681 


2.248 
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0.847 


0.912 


0.850 


0.911 


7 


3.581 


3.674 


11.64 


11.56 



dot-dashed lines for SLy and FPS, respectively). It illus- 
trates the accuracy of the fit l|14|) . 

Now, n(p) can be easily obtained from Eq. (|3J). Doing 
this and substituting P(p') in the integrand of Eq. J3J 
from Eq. (|14|l . we recover the original tabular values with 
maximum difference < 0.4% for FPS and < 0.12% for SLy. 

In some applications, it may be convenient to use n as 
independent variable, and treat p and P as functions of n. 
For this purpose one can use the following fit: 



nmo 



pin P2 + p^n P4 

(1 +P5Tl) 2 

n 



8 x 10- 6 + 2.1 n - 585 



/ (-P6(logn + p 7 )) 

/ (p 6 (logn + p 7 )) , (15) 
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Fig. 3. Dependence of p on n. Upper panel: rarefied data 
(symbols) and the fit l|15|l (lines); lower panel: relative 
difference between the data and fit. Filled dots and solid 
line: SLy; open circles and dot-dashed line: FPS. 



where n is in fm 3 . The inverse fit n(p) is given by the 
formula 

X q x x q2 + q 3 x q4 

- = 1 + —rr l s3-/o(56(97-iogp)) 

n (l + q 5 xy 



x 



8 x f0- 6 + 2.1a; - 585 



/oOetlogp- q 7 )), (16) 



where x = p/mo and p is in g cm' 3 . Coefficients pi and 
q. t of the fits (JTHJ and (JTSJl are given in Table [3 The 
fractional fit residual in p and n varies from ^ 10 -9 at 
p 10 g cm -3 to a fraction of percent near the maxi- 
mum p ~ 10 16 g cm -3 , while the difference (p — nmo) is 
approximated with the typical accuracy of a few percent. 
In Fig. |3 (upper panel) we plot the ratio p/(nmo) against 
n for the two EOSs. We see that in the crust, at n ^ 0.1 
fm -3 , p is to a good accuracy proportional to n. The lower 
panel, which shows the relative difference between fitted 
and tabulated values of p(n), confirms the accuracy of 
Eq. (JTSJl. 

It should be stressed that thermodynamics requires 
Eq. to be satisfied exactly. To achieve this, one should 
not totally rely on the fits ljl5[) and Ijlfijl : otherwise ther- 
modynamic consistency will be violated on the scale of 
these fits' errors (a fraction of percent). So, if p is used as 
an input, then rtb(p) should be calculated from Eq. 
Alternatively, if the input is Tib, then, having calculated 
Pfit(nb) from Eq. (THJ) and P(n h ) = P(pat(n h )) from 
Eq. (|14J) . one should refine p(rib) using relation (J2J). 



Fig. 4. SLy and FPS EOS for rotational configurations. 
Upper panel: rarefied data to be fitted, calculated accord- 
ing to Eq. I|14[) (symbols) and the fit l|17|l (lines); lower 
panel: relative difference between the data and fit. Filled 
dots and solid line: SLy; open circles and dot-dashed line: 
FPS. 

3.2. Rotating stars 

As explained in Sect. 12.21 for rotating stars it is most useful 
to parametrize density and pressure as functions of the 
pseudo-enthalpy H, which can be written in terms of the 
enthalpy per baryon h. Let us define r\ = h/moc 2 — 1. 
In view of relation ©, the function which we intend to 
parametrize, £(7/), is not independent of the function £(£) 
parametrized by Eq. (|14fl . In order to fulfill Eq. JSJ) as 
closely as possible, we first calculate ri(£) using Eqs. I|14f> 
and ©, and then find the inverse fit £(77). The best fit 
reads: 



bi + 6 2 log r\ - 



M b4 



,/o(Mlog77- M) 



1 + M. 

b 8 + b 9 log 77 + (b w + frn log y) {b\2Tl) 7 
1 + &1377 + {bi2V) 7 
x/o (K(b 7 - log 77)) + b 14 f (b 15 (b 16 - log 77)) 



(17) 



where the parameters bi are given in Table The com- 
parison of the fit and the data is presented in Fig. The 
typical fit error of p, according to Eq. fTTjl . is about 1% 
at 77 <; 10 (corresponding to £ ^ 3), and the maximum 
error < 4% occurs near the neutron drip and crust-core 
phase transitions. 

When used in combination, the fits l|14(l and (|17|) . to- 
gether with Eq. @ or Eq. (|16|l give the parametrizations 
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Table 3. Parameters of the fit ([TT| 
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bi (FPS) 


bi(SLy) 
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6i(FPS) 


6* (SLy) 
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5.926 


5.926 
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11.97 


34.96 
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0.4704 


0.4704 


10 


15.432 


15.328 
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19.92 


20.13 


11 


0.6731 


0.621 
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0.2333 


0.2347 


12 


49.4 


63.1 
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2.63 


3.07 


13 


11.47 


68.5 
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54.7 


97.8 


14 


1.425 


2.518 
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-1.926 


-2.012 


15 


3.0 


2.6 
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36.89 


89.85 


16 


0.913 


1.363 



of p(H), P(H), and n(H) needed for calculations of the 
stationary rotating configurations. In this case, the func- 
tion P(H) = P(p(H)) obtained using Eqs. JJJt and EM , 
reproduces the tabular values with a typical discrepancy 
of 1-2%, with a maximum within 10% near the crust-core 
boundary. 

The remark on the thermodynamic consistency, made 
at the end of Sect. 13. ll applies also here: one should refine 
either or p fitted values, using the exact relations Q 
or 10). 

4. Adiabatic index 

An important dimcnsionlcss parameter characterizing the 
stiffness of the EOS at given density is the adiabatic index, 
defined by 



n dP 
P dn 



1 



P 



p dP 
pc 2 P dp 



(18) 



Using our fit 1|14[). we obtain the analytical expression 
a 2 — aia 4 + 3a 3 £ 2 + 2a 3 a 4 £ 3 



p dP _ d( 
P dp ~~ d£ 



(l + a 4 £) 5 



-0-5 7- 7 fo{a 5 {a 6 - Q) 

1 + a 4 £ 

4 

■ jo( a 4j+l ( a 4i+2 — 0) [ a 4i 
i=2 

+a 4i+ i(a 4i _i + a 4i £)/ (a 4l+1 (£ - a 4i+2 ))] 



(19) 



Different regions of neutron-star interior are character- 
ized by distinct behavior of T, displayed in F ig. El Precise 
values of T calculated by iDouchin fc Haensell l|2f)f)lh are 
shown by dotted line, and the analytical approximations 
according to Eqs. I|14f> . (|19l) by solid line (for SLy) and 
dot-dashed line (for FPS). 

In the outer crust, the value of T depends quite weakly 
on density. At 10 8 g cm -3 Ss p & 3 x 10 11 g cm -3 , this 
value would be T w 4/3 if A, Z values were fixed, because 
in this case P is mainly determined by the pressure of ul- 
trarelativistic electron gas which behaves as oc (Zp/A) 4 / 3 . 
For example, T « 4/3 within each shell with con s tant A 
and Z of the EOS derived bv lHaensel & Pichonl l(l994h . 
However, the compressible liquid drop model used by 




Fig. 5. Adiabatic index T for SLy and FPS EOSs. Solid 
line: analytical approximation (SLy); dotted line: precise 
values (SLy); dot-dashed line: analytical approximation 
(FPS). 



IDouchin fc Haensell l|200lj) effectively smoothes the dis- 
continuities caused by transitions from one to another 
(A, Z)-species with increasing density, which leads to an 
effective continuous increase of the A/Z ratio and corre- 
sponding decrease of T, seen in Fig. [3] 

A dramatic drop in T occurs at neutron drip thresh- 
old, which corresponds to strong softening of the EOS. 
The analytical expression i|19|) somewhat smoothes this 
drop. The behavior of T in the inner crust results from 
an interplay of several factors, with stiffening due to in- 
teraction between dripped neutrons, a softening effect of 
neutron-gas - nuclear-matter coexistence, and the soften- 
ing Coulomb contribution. 

At the crust-core interface, matter strongly stiffens, 
and r jumps from w 1.7 to w 2.2, which results from 
disappearance of nuclei. The analytical approximation 
smoothes this jump also, though reflects the stiffening. 
This approximation is also smooth across a small dis- 
continuous drop of r at p w 2 x 10 14 g cm~ 3 , where 
muons start to replace a part of ultrarelativistic electrons. 
However, the electrons and muons give only minor contri- 
bution to the pressure (and therefore behavior of T) in the 
core, because the main contribution comes from interac- 
tions between nucleons. 

5. Summary, discussion, and conclusion 

Analytical representations of the EOSs used in the modern 
2-D and 3-D simulations of neutron star dynamics have 
many important advantages over the tabulated EOSs. 
Analytical EOSs do not require any interpolations, are 
thermodynamically consistent, and allow for a very high 
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precision of calculations. In the present paper we con- 
structed analytical representations, in terms of the contin- 
uous and differentiablc functions of a single chosen vari- 
able, of the SLy and FPS EOSs. For these analytical rep- 
resentations, the thermodynamic relations are exactly sat- 
isfied at any point. Two choices of the independent vari- 
able were considered. The first one is p; function P(p) 
given by Eq. (|14fl fits the original tables in the density 
interval 10 5 g cm -3 < p < 10 16 g cm -3 within typical er- 
ror of 1-2%. Function n(p) can be calculated either from 
Eq. (PJ) to satisfy exactly the first law of thermodynam- 
ics, or from the fit (|16[) with a typical error ~ 0.1%. A 
variant which ensures the same accuracy is to choose n 
as an independent variable and calculate p(n) from the fit 
lfl5|l and P(p) from Eq. (11} . The other choice of the in- 
dependent variable is to use the pseudo-enthalpy H . This 
choice is particularly advantageous for applications to 2-D 
and 3-D numerical simulations of neutron star dynamics, 
such as rotation and inspiraling stage of the evolution of 
relativistic neutron-star - neutron-star binary. We repre- 
sented both EOSs by the continuous and differentiable 
functions P(H), p(H), and n(H), where p(H) is given by 
Eq. (|17fl with typical accuracy within a few percent, while 
P(H) and n(H) are calculated from the functions P(p) 
and n(p), respectively. Differentiation of P{p) then yields 
analytical representations of the adiabatic index for the 
SLy and FPS EOSs; this quantity is important, for exam- 
ple, for numerical simulations of dynamics of the neutron- 
star - neutron-star system at the inspiral phase. 

The quality of our analytical representations of the 
EOSs was tested by evaluation of the virial identi- 
ties G RV2 llBonazzo lai Il973t iBonazzola & Courgoulhorl 
Il994|) and GRV3 ((Gourgoulhon fc Bonazzolalll994l) in the 
numerical simulations of the 2-D stationary rotation of 
neutron stars. GRV2 and GRV3 are integral identities 
which must be satisfied by a stationary solution of the 
Einstein equations and which are not impo sed in the nu- 
merical procedure (see iNozawa et al.lll998l for the details 
of computation of GRV2 and GRV3). For rotating con- 
figurations we get GRV2, GRV3 - 1(T 6 - 1(T 5 , which is 
excellent. 

The subroutines for the numerical applica- 
tions of analytical EOSs (Fortran and C++ ver- 
sions) can be downloaded from the public domain 
http : //www . iof f e . ru/astro/NSG/NSEDS/. 
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